Nonequilibrium roughening transition in a simple model of fungal growth in 1 + 1 
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We introduce a simple model of yeast-like growth of fungi colonies, which exhibits a continuous 
roughening transition far from equilibrium from a smooth (a — 0) to rough phase (a = 1/2) in 1 + 1 
dimensions. At the transition some scaling exponents are calculated by mapping the problem onto 
a directed percolation process. The model reproduces the roughening transition observed in some 
experiments of fungal growth. 
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Much effort has been devoted to the search for the ba- 
sic principles governing the pattern formation in living 
organisms. Among all the phenomena of formation and 
development of complex structures involving living or- 
ganisms, the growth of bacteria and fungi colonies has 
attracted a considerable amount of work in recent years 

Besides the formation of patterns, which is associ- 
ated with the existence of some unstable modes, spatio- 
temporal scale invariance is also commonly observed de- 
pending on the environmental conditions ||7|-|l0|. In the 
particular case of fungi colonies, the morphology may be 
well classified in hyphal and yeast-like growth The 
former corresponds to multicellular growth and fractal 
filamentous patterns form due to the existence of highly 
cooperative behavior of the individual cells. However, 
yeast-like growth occurs in solidified media and in this 
case the colony is a very compact object. The front 
of the colony usually becomes rough (i.e. the interface 
width diverges with the linear size of the system) and 
its dynamics is completely characterized by the critical 
exponents, which describe the scaling properties of the 
interface fluctuations |14|. 

In recent experiments [fl2[ with the yeast Pichia mem- 
branaefaciens on solidified agarose film, several morpho- 
logical transitions in 1 + 1 dimensions have been reported. 
In these experiments different front morphologies were 
obtained depending on the concentration of polluting 
metabolites. Also transitions from rough to flat inter- 
faces were observed fll5| , although never studied system- 
aticly. In the case of Bacillus subtilis, the existence of 
completely smooth and Eden-like morphologies (among 
others) depends on the agar hardness and nutrients con- 
centration pT| . Our attention in this Letter is focused 
on this roughening transition. 

In 1 + 1 dimensions a surface under equilibrium condi- 
tions is rough at any finite temperature , because ther- 
mal fluctuations always make the flat surface entropically 
unfavorable. Only in higher dimensions can equilibrium 
interfaces undergo a phase transition from a rough to a 
smooth interface at some critical temperature. The situ- 
ation in nonequilibrium is much richer because detailed 



balance is not required and growth processes may ex- 
hibit a roughening transition even in 1 + 1 dimensions, 
although examples are few |]l6|-|l8| . More generally, phase 
transitions in nonequilibrium 1 + 1 dimensional systems 
have usually been observed in systems with absorbing 
states. Thus, it is of great interest to find models far 
from equilibrium which do not possess absorbing states 
but still display a phase transition. 

In this Letter we introduce a new class of models di- 
rectly inspired by yeast-like fungal growth in 1 + 1 dimen- 
sions. The model has no absorbing states and exhibits a 
nonequilibrium roughening transition between a flat and 
a rough phase. Some scaling properties of the model at 
the transition can be related to directed percolation. We 
obtain that close to the transition a diverging correla- 
tion length £ ~ |e| _l/ appears, where e is the distance to 
the critical point, and the critical exponent is in good 
agreement with directed percolation in 1 + 1 dimensions, 
v ss v±_ ss 1.10. We find that, in the rough phase, the 
model belongs to the Edwards- Wilkinson [[19| universal- 
ity class. 

Model. - We study a class of models motivated by the 
growth of colonial organisms like fungi and bacteria. The 
model is defined on a 1 + 1 dimensional lattice in which 
every site can have two different states: occupied or va- 
cant. Growth of the colony occurs because of the division 
of individual cells, so only nearest neighbors of occupied 
sites can potentially become occupied. The basic idea 
now is that cell division is less likely in young cells, so 
the probability for a vacant site i of being occupied in the 
next time step has to increase with the total age A4 (t) of 
the occupied nearest neighboring sites of that vacant site 
i. At every time step a growth probability 



Pi(t) = F(6 m Mt)) 



(1) 



is assigned to all vacant sites that are nearest neighbors 
of occupied sites (note that Ai — at the remaining va- 
cant sites). F(x) is any monotonous increasing function 
that satisfies F(0) = and F(oo) = 1. For defmiteness, 
in our simulations we simply used F(x) = tanh(ir). 9 
is an external parameter which controls the growth rate. 
The variation in the aging times among the individual 
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cells in the colony is modelled by the parameter r\i. This 
random variable is quenched, uncorrelated and uniformly 
distributed in [1 — A, 1 + A]. The simulations are per- 
formed on a triangular lattice with open boundary con- 
ditions, parallelly updating all the sites according to the 
growth probability (Q) . The model produces a very com- 
pact bulk and irrelevant overhangs at the scale of the 
lattice spacing may appear on the growth front depend- 
ing on the tunning parameter 9. 

From the point of view of the biological process, this 
model has to be considered as a simplified version of a 
more realistic description. The experiments mentioned 
above |L2| have shown that yeast-like growth of fungi 
strongly depends on the concentration of nutrients and 
inhibitors in the medium. In a more complicated model 
p0| the effect of inhibitors can be studied by including 
an inhibitors diffusive field Ci(t) which is coupled with 
the growth probability at every site. The probability for 
a vacant site i of becoming occupied is then given by 
Pi{t) — F(9r]iAi(t)exp(—Ci(t)/co)), where cq is a thresh- 
old concentration above which the effect of inhibitors be- 
gins to be relevant. The model we are presenting here 
corresponds to the limit cq — > oo. 

Continuous Roughening Transition.- Despite its sim- 
plicity the model has a non trivial behavior. For small 
values of the parameter 9, the fungus front generated by 
the model becomes rough and its dynamics is character- 
ized as usual by two kinetic roughening exponents [fl4|| . 
In contrast, for large values of 9, the interface becomes 
flat even when started from a non smooth initial condi- 
tion. This phase transition occurs for any value of the 
disorder intensity A. We found that these two phases 
are separated by a critical threshold which depends on 
the degree of disorder and is at 9 C = 0.183 ± 0.003 for 
the case of no disorder A = 0. In the following we re- 
strict ourselves to the homogeneous case A = although 
similar results are found for any < A < 1. 

Let us first study the rough phase, 9 <C 9 Cl in some 
detail. Denoting by h(x, t) the front height at time t 
and substrate position x, the interface width in a system 
of size L, W(L,t) = {[h(x,t) - Ht)} 2 ) 1 / 2 (where overbar 
and (• • •} denote spatial average and over realizations re- 
spectively), measures the magnitude of the interface fluc- 
tuations. The width is expected to satisfy the dynamic 
scaling ansatz 21 13]: 



W(L,t) 



t 13 if t < L z 
L a if t > L z 



(2) 



where a is the roughness exponent and z = a/ f3. 

We have performed simulations in systems of differ- 
ent sizes to determine the universality class of the model 
in the rough phase. Figure 1 shows that sufficiently far 
from the transition (6 — 0.01) the scaling behavior of the 
interface is given by Eq.(g). The kinetic roughening ex- 
ponents a = 0.46 ±0.05 and j3 = 0.24 ±0.02 indicate that 



the model belongs to the Edwards- Wilkinson universality 
class fy§. 

Approaching the roughening transition from the rough 
phase a new diverging length £ ~ e~ u , where e = 9 C — 9, 
appears and enters the scaling of the width. Simula- 
tions close to the transition (e — > + ) were performed in 
systems of different sizes. In Fig. 2 results for W(L,t) 
vs. time corresponding to L = 900 and average over 
50 realizations are plotted. We found that there ex- 
ist two different regimes separated by a crossover time 
t c ~ e -7 , which diverges as e — + + . Figure 2 shows 
that , near and below the transition, W(L,t) scales as 
e' KF for times t < t c < L z and W(L,t) ~ t !3 e KR for 
t c <C t <C L z . kf (kjj) stands for the scaling exponent in 
the flat (rough) regime of the surface width as the tran- 
sition is approached. For times t~^>L z the interface gets 
to saturation as in ([I). 

The scaling behavior in the intermediate time regime 
before saturation suggests the following scaling ansatz for 
times t <C L z : 



W(L,t,e) 



r g(t/t c 



(3) 



where the scaling function is g{u) ~ const if u -C 1 and 
g(u) ~ u 13 for u 1. The scaling relation 



(4) 



between critical exponents must also be fulfilled. The in- 
set of Fig. 2 shows the data collapse for the exponents 
7 = 1.73 ± 0.03 and up — 0.02 is consistent with Eq.(|^) 
over a range of more than four decades. Due to the very 
small value of Kp a logaritmic dependence W ~ log(e) 
cannot be ruled out. Using the scaling relation (^) we 
find kr ~ 0.41, which is to be compared with the the 
measured value kr ~ 0.42 ± 0.01. 

As seen in some other 1 + 1 dimensional growth models 
exhibiting a roughening transition [|l7|jl8| , some scaling 
properties are related to directed percolation (DP) [E2|. 
Our model turns out to be a model with a maximal veloc- 
ity. Small values of 9 produce small growth rates and the 
front propagates with a finite velocity. If the parameter 
9 is increased, the probability for a site of growing (i.e. 
becoming occupied) is also raised as follows from Eq. ([[]) . 
A given site reaches the maximal velocity when its proba- 
bility of growing at the next time step is one and then no 
further increase of the propagation velocity at that site 
can occur [p3| . The transition from rough to flat surface 
is a DP transition. The flat phase corresponds to the 
active DP phase whereas the rough phase corresponds 
to the non-active DP phase. The sites at the highest 
level h — v max t are active sites of a DP process. Below 
the transition (rough phase) there are no sites growing 
at maximal velocity and the surface mean velocity is less 
than v m ax ■ Above the transition (flat phase) , on the con- 
trary, there is a finite density of sites advancing at v max . 
At 9 = 9 C the sites which grow with maximum velocity 
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constitute a percolation cluster directed along the growth 
direction, as depicted in Fig. 3. 

The mapping to DP allows the derivation of several 
scaling exponents. The correlation length £ is identified 
with the correlation length of the DP clusters which 
is transversal to the growth direction of the interface, and 
v = u± w 1.10 in 1 + 1 dimensions. This implies that the 
crossover time t c ~ £ z in Eq. (||) is given by t c ~ using 
the fact that z' — V\\/v± for DP. Our estimation of the 
crossover time exponent 7 = 1.73 ± 0.03 (see Fig. 2) is 
in excellent agreement with the DP exponent vu w 1.73. 

Order Parameters.- We now define a convenient order 
parameter for the rough phase (9 < 9 C ) as V = v max — v 
with v = dth(t) the mean front velocity. The variable V 
measures the fraction of the surface points propagating 
at velocity lower than the maximal one. V is different 
from zero only in the rough phase and its scaling be- 
havior can also be related to properties of the associated 
DP process. In the rough phase V is given by the in- 
verse of the typical time that the DP correlations survive 
t ~ £|| , which is the only characteristic time scale. So, 
V ~ 1/t ~ e^n close to the transition. In our simulations 
we measured V ~ e 1 - 74e± 0- 018 j as shown in Fig. 4, in 
excellent agreement with the simple DP scaling picture. 

As for the flat phase (9 > 9 C ) the density of sites at the 
maximal height is an appropriate order parameter [ jl7| , 
which corresponds to the density of active sites in DP. 
However, more interesting here may be the use of ap- 
propriate order parameters which take into account the 
symmetries of the model in order to study the possibil- 
ity of spontaneous symmetry breaking (SSB) in the flat 
phase. It has been shown that SSB may take place in 
nonequilibrium situations under certain conditions even 
in 1 + 1 dimensions p|p4| . 

The dynamical rules of our model are invariant under 
translation of an integer number of levels in the growth 
direction. The dynamics do not favor odd nor even 
heights. Thus it is adequate to define a magnetisation- 
like order parameter 

M(t) = ^(-l) h ^\ (5) 
»=i 

which is no conserved by the dynamics of the model. This 
order parameter can be used to quantify the symmetry 
breaking which may take place in the flat phase |M . 

In our simulations we found that in the rough phase 
(M) = and also (\M\) = 0, where (• • •} denotes an 
average over realizations. In the rough phase the inter- 
face explores many height levels and the contribution to 
the magnetisation at different sites average out. On the 
contrary in the flat phase (M) = but (\M\) 7^ due 
to the fact that most heights are at the highest level. 
However, SSB was not observed in our model. For ev- 
ery realization M(t) gets to a stationary state consisting 
of flips between + |M S | and — \M S \ at almost every time 



step. In Ref. [JISJ] a growth model with sequential up- 
dates has recently been studied in which SSB was found. 
In that model, SSB in the flat phase is related to the fact 
that the interface eventually becomes pinned at a certain 
height and the velocity is actually zero in the thermo- 
dynamic limit. Only in a finite system will the interface 
have a finite velocity which vanishes exponentially with 
the system size. This ensures the breaking of the symme- 
try between odd and even heights in the thermodynamic 
limit. The final interface will be pinned at an odd or 
even height, depending on the initial conditions, break- 
ing the symmetry of the dynamics. On the contrary in 
our model the velocity v max in the flat phase is always 
finite (even in the thermodynamic limit L — > 00) not al- 
lowing SSB. As occurs in the model reported in Ref. Jl?]] 
the smooth phase does not survive in our model under 
sequential dynamics, due to the random walk nature of 
sequential updates. 

The behavior of the order parameter of the flat phase 
near the transition gives a new scaling exponent 77, 
which describes its decay as the transition is approached 
from above (\M\) ~ (— e) n . Our simulations indicate that 
77 ~ 0.50 which can be compared with 77 = 0.55 ± 0.05 
previously measured for a different model (ll| . 

Conclusions.- We have studied a model for fungal 
growth in which the likelihood of growth at a given site 
depends on the local environment. Vacant sites sur- 
rounded by older cells are more likely to become occu- 
pied than sites neighboring younger cells. The model 
exhibits a roughening transition in 1 + 1 dimensions as 
the growth rate is increased. We have found that the 
model belongs to the Edwards- Wilkinson universality in 
the rough phase. At the transition some critical expo- 
nents have been calculated by mapping the problem onto 
a DP process. Transition from a completely flat to rough 
interface has been observed by T. Sams et. al in experi- 
ments [p"5[ with colonies of the yeast Pichia membranae- 
faciens on agarose film and by Wakita et. al. p"T) in the 
bacterium Bacillus subtilis. 
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FIG. 1. Interface width of the front vs. time in 
the rough phase (6* = 0.01) for different system sizes 
L = 16, 25, 32, 50, 64, 75, 100 . The slope of the dotted line 
corresponds to the time exponent j3 = 0.24±0.02. In the inset 
the values of the width at saturation are plotted vs. system 
size. The dotted line fits the data and gives the roughness 
exponent a = 0.46 ± 0.05. 



FIG. 2. Time behavior of the width in the rough phase 
near the transition. Curves correspond to different distances 
to the critical point, form e = 0.083 (top) to e = 0.003 (bot- 
tom). Data are collapsed in the inset following Eq.(3) for 
7 = 1.73 and kf = 0.02, where the dotted line has slope 0.25. 




% = 1 .0029 o 

FIG. 3. Clusters of sites growing at maximal velocity. Ev- 
ery pixel represents the position (x, h) of a site which moved 
at Umax and configurations at intervals of 5 time steps are 
shown. The system size is L = 500 and time advances up- 
wards. 
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FIG. 4. Scaling of the order parameter V in the rough 
phase near the transition. The line is a least-square fit of 
the data and has a slope 1.746 ± 0.018 in good agreement 
with vu . 
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